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A RISK ANALYSIS FOR A SYSTEM STABILIZED BY A 
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Abstract. We formulate and analyze a multi-agent model for the evolution 
of individual and systemic risk in which the local agents interact with each 
other through a central agent who, in turn, is influenced by the mean field of 
the local agents. The central agent is stabilized by a bistable potential, the 
only stabilizing force in the system. The local agents derive their stability 
only from the central agent. In the mean field limit of a large number of 
local agents we show that the systemic risk decreases when the strength of the 
interaction of the local agents with the central agent increases. This means 
that the probability of transition from one of the two stable quasi-equilibria to 
the other one decreases. We also show that the systemic risk increases when 
the strength of the interaction of the central agent with the mean field of the 
local agents increases. Following the financial interpretation of such models 
and their behavior given in our previous paper (Gamier, Papanicolaou and 
Yang, SIAM J. Fin. Math. 4, 2013, 151-184), we may interpret the results of 
this paper in the following way. From the point of view of systemic risk, and 
while keeping the perceived risk of the local agents approximately constant, 
it is better to strengthen the interaction of the local agents with the central 
agent than the other way around. 

Mean Field Models, Dynamic Phase Transitions, Systemic Risk 

1. Introduction 

In recent years, interacting particle systems have been extensively used to model 
financial systemic risk for complex, inter-connected systems. An interacting particle 
system with binary risk variables is considered in [3] and the law of large numbers, 
central limit theorem and large deviation principle are derived for this model. An 
interacting particle system of diffusion processes is used in [3] to model the interbank 
lending system. In [3], a model simplified from the one in [9] is considered, in which 
each agent can control the lending flow rate and optimizes the individual objective 
function, and thus the system can be put in the framework of mean field games. In 
HD, the authors use interacting Bessel-like diffusion processes to model systemic 
risk and establish a large deviation principle. In Emm, we consider an interacting 
particle system with a bistable potential and we use the large deviation principle 
to explain that the overall systemic risk may increase while individual risks are 
decreased. The large deviation principle in TT) TJ| is solved numerically in [171 . In 
PP, the authors consider interacting jump-diffusion processes modeling interbank 
lending and borrowing and prove the weak law of large numbers (LLN) of the 
empirical measure as the number of individuals goes to infinity, and define systemic 
indicators based on the LLN result. In EH Em OH El], the authors model large 
portfolios and default clustering and derive the law of large numbers, fluctuation 
analysis and large deviations. 
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In our previous work [10], we used an interacting agent-based, mean-field model 
to show that individual risk may not affect systemic risk in an obvious way. That is, 
each agent may have relatively low individual risk by diversification through risk¬ 
sharing while the overall, systemic risk is increased as a result of diversification. We 
considered the following model that was studied extensively before by [51 [61 EH [7] : 

(1) dxj(t) = —hV'(xj{t))dt — d(xj(t) — XN{t))dt + ad,W {, j = 1,... ,7V, 


where Xj ( t ) represents a risk variable for agent j at time t and N is the number of 
agents. The potential V{x) = \x 4 — \x 2 is taken to be bistable with two stable 
states ±1, and the constant h > 0 quantifies intrinsic stability for each agent. We 
define —1 as the normal state of an agent and +1 as the failed state. The empirical 
mean Xjv(t) ■= Y^jLi x j(t) is the mean risk, and the constant 6 is positive so 
that Xj tends to stay close to xjy. The standard Brownian motions are 

independent and model external risk factors, with er > 0 their strength. 

It was shown in [5] that the empirical measure UN(t,dx) := (t) (dx) 

converges weakly in probability to u(t,dx) = u{t 1 x)dx , the weak solution of the 
nonlinear Fokker-Planck equation: 


|«=4r<xH-4 


yu{t,dy) - x 


1 2 d 2 
2 ct dtf u ' 


starting from u(0, dx) = lirn/v^oo Un(0 , dx) (provided the weak limit exists). Given 
h and 9 , for sufficiently small a , u{t , x) has two equilibria u±^ b (x) := lim^oo u(t, x), 
where x^v (0 converges to either £5 > 0 or — ^ as t —> 00 , depending on the initial 
condition. Thus we define u e _^ b as the normal state of the system and u [as the 
failed state of the system. 

Given that N is large but finite, and U]y(0,dx) « u e _^ b {x)dx, we showed [TOJ 
Theorem 6.2 and Corollary 6.4] that by using the large deviation principle in [ 6 j 
and assuming that h is small, the systemic risk, defined as the probability of the 
transition of Uisr(t,dx) from u e _^ b {x)dx at time 0 to u e ^_^ b (x)dx at some time t < 
T < 00 has the following exponentially small but nonzero value: 

( 2 ) 

P (Un( 0, dx) ~ u e _£ b {x)dx, UN{t,dx) ~ u e + ^ b {x)dx t <T < 00 ) exp f 


where 


= 



6 (a 2 \ 1 - 2(a 2 /20) 

1 - 3 ( ct 2 / 26 ») 


1 + h ^ [ 29 


0{h 2 ). 


Fluctuation analysis on 0 m Lemma 6.5], shows that the risk of each agent has 
the form Xj(t ) = — 1 + Zj(t) and lim^oo Var^j-(t) < | -s- Thus, the quantity can 
be considered as the individual risk for each agent. 

We then see that if the strength of the external risk a 2 is increased, either 
because the agents are more risk-prone or because the economic environment is 
more uncertain, then the agents can increase 9 , the risk-diversification parameter, 
so that that their individual risk is still low. However, from the analysis of the 
systemic risk 0 we see that the systemic risk is increased when er 2 increases even 
if the individual risk er 2 / (29) is very low: there is a systemic level effect of a 2 that 
cannot be observed by the agents and it tends to destabilize the system. 
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In this paper, we extend the previous model 0 by introducing a central agent 
with the risk variable Xg N \t). The model we study in this paper is given by 

1 N 

(3) dx^ = —hoVo(xQ N ^)dt - d 0 (x < ' 0 N) - x N )dt + -jLdW?, x N = — x j: 

A j —1 

(4) dxj = — hV'{xj)dt — d[xj — Xq N ^)dt + adW ?, j = 1,..., N. 

Here Vo(x) and V(x) are potentials with two stable states and in this paper we 
again assume that Vg(x) = V(x) = \x A — \x 2 with the stable states ±1. The 
parameters h 0 , h > 0 are the strengths of intrinsic stability of the central and local 
agents, respectively. The parameters 9q,9 > 0 determine the strength of the mean- 
field interactions. The central agent x [^ is intrinsically stable when ho > 0 and 
may be destabilized through a mean field interaction with the local agents where 
9 0 > 0. Depending on whether h > 0 or h = 0, the local agents {aL,}^ are or are 
not intrinsically stable. They may be stabilized through their interaction with the 
central agent Xq N \ The independent, standard Brownian motions {IF/I^Lq model 
the external risk for the central and local agents. We note that the normalization 
factor 1 /y/N in |s| makes Xq N ' 1 and xn have external risks of comparable size for 
N large, and we will assume that <to < cr or <to = 0 since we want the central agent 
to operate with less risk than the local agents. 

In the regime of no cooperation, 9g = 9 = 0, the central agent and the local 
agents are independent of each other and Kramers’ large deviation law states that 
when (To and a are small, the probabilities of transition from one stable state to 
the other within the time interval [0,T] are proportional to Texp(— 2hoVo(t)) / cr'g) 
and T exp(— 2hV(0)/a 2 ), for the central and local agents, respectively. We want to 
analyze stabilization effects in the cooperative regime 9g,9 > 0. 

In this paper, we will assume that the intrinsic stability of the local agents, 
h, is exactly zero, while we only assume that h is small in m ■ Because of this 
simplifying assumption, instead of considering the pair (xg N \t), ^Xj(t){dx)) 

as a scalar and a measure-valued process, we can simply consider (xg N \t), Xjv(t)) 
as a two-dimensional process and get results that are more detailed than it was 
possible in the setup of HD]. First, we compute numerically the minimizing path 
for the associated large deviation problem, and we are able to explore how the 
various parameters affect the agents’ fluctuations and the systemic risk. We also 
recover the main result in [lOj . that is, that the systemic risk is increased, with the 
local risks kept fixed, if we increase er 2 and 9 with the ratio a 2 /9 fixed. Another 
result is that because we assume that 0 = h. < h 0 and a 0 < er, the central agent is 
more stable than the empirical mean of the local agents. In this setting, we find 
that 9q and 9 tend to play opposite roles: higher 9q increases the systemic risk as 
we force the stable term Xq to be close to the relatively unstable term x, but 
on the other hand, increasing 9 lowers the systemic risk as x tends to be close to 
Xq N \ This is the main result of this paper. The third result here, for a case not 
considered in the previous paper, concerns the introduction of optimal controls for 
the local agents. We use optimal control theory and find that the use of controls 
amounts to replacing 9 by an effective one that is larger, and thus it reduces the 
systemic risk. 
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This paper is organized as follows. In Section [2] we state the mean field limit of 
the pair (xq N \t), $xj(t)(dx)) as N —» oo. We then discuss the equilibria of 

the limit Fokker-Planck equation. In Section [3] we analyze the special case where 
h is exactly zero. In this case, explicit solutions of the fluctuation analysis can be 
obtained, and we have a large deviations principle for (xq N \t),XN(t)) using the 
Freidlin-Wentzell theory. In Section [4] we give the formal large deviation principle 
for the empirical measure (xq^ (i), Y^J=i d Xj (t)(dx)) that is necessary when ft > 0 . 
We do not use this general formulation but we do show that the large deviation 
problems for (xg N \t), Xjv(t)) and (xg N \t), J2‘j=i d Xj (t){dx)) are the same if ft = 0. 
In Section [5] we formulate a control problem for the local agents in Q and use 
optimal control theory to analyze the effect of the control on the system. Finally, 
in Section [ 6 ] we present results of extensive numerical simulations. The technical 
details of the proofs are in the appendices. 


2. The mean field limit of a large number of local agents 

We begin by recalling the main results of mean field limit theory as they apply 
to problem ©,@, in the next section, and then discuss the equilibrium solutions 
of the limit, non-linear Fokker-Planck equation. 


2.1. The non-linear Fokker-Planck equation. The stochastic model M is 
a simple extension of the model in [5] [12] (see also (23U22UI8UI6]). We let Mi (M) 
denote the space of probability measures endowed with the metric of the weak 
convergence, and C([0, T], Mi(R)) the space of continuous Mi (M)-valued processes 
in the time interval [0,T] endowed with the maximum distance in [0,T]. In the 
limit N — > oo, the pair (xg N \t), &Xj{t){dx)) converges in (R,Mi(M)) to 

(yo(t),p(t,x)dx) in probability, the weak solution of the nonlinear Fokker-Planck 
equation and ordinary differential equation 

(5) = -h 0 V'(yo) - 9o(yo - J xp(t, x)dx s j , 

(6) x) = h-^- [V\x)p(t, x)] + 9-^\(x- yo(t))p(t, x)] + x), 

with the initial condition z/ 0 (0) = linqv-ioo Xo^(0) andp(0, dx) = lim^-ioo SjLi d x 
given that the limits exist. Equivalently, we can characterize the pair (yg(t),p(t, x)dx) 
by noting that p(f,x) is the transition probability density of the process X t , the 
solution of 

^yo = - h 0 Vo(y 0 ) - 9 0 (yo ~ EX t ), 

dX t = —hV'(X t )dt — 9{X t — yo)dt + (rdW t , 


where Wt is a standard Brownian motion. In addition, if ft = 0 and y(t) := E(A' t ), 
then ( y 0 (t),y(t )) satisfies 


(7) 

( 8 ) 


^yo = -h 0 Vg{y 0 ) - 9 0 (y 0 - y), 
4y = ~9(y- yo). 


i(o )(dx), 
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2.2. Equilibrium states. Given the existence of a stationary state (yg,p e (x; yg)) := 
(lim t 

—>•00 i/o(i),lim t _ i . 00 p(t,a;)), it satisfies 


(9) 


p e (x; 2/0) = 


z(y e 0 ) 


exp 


2 hV ( x ) + 0(x — yg)' 


which is obtained from ([ 6 ]), and satisfies the consistency equation 

ho, 


( 10 ) 


xp e (x\y e 0 )dx = Vo + ^-^'(yg), 


obtained from (§. If h = 0, then p e (x ; yg) is a Gaussian density function, given by 
& with mean yg and ( fl0| ) implies Eg(yg) = 0. Therefore yg = ±1. The equilibrium 
states for the system are determined by the equilibrium states of the central agent. 
Indeed, if the central agent takes the equilibrium value yg = —1, then the individual 
agents take a Gaussian distribution with mean —1 and variance a 2 /(20): 


( 11 ) 


p e (x) = 


exp - 


0( x + 1 )' 


When h is positive but small, we let yg° = ±1 and therefore Vg(yg°) = 0 with 
Vg"(yg°) > 0. It is then possible to find an equilibrium state yg - , resp. yg + , close 
to yg° = — 1 , resp. yg° = 1 , and we have 

Vo = Vo° + h y o 1 + o(h), 


with 

el _ 0 o fe- ex2 /° 2 V'(y% 0 +x)dx 

Vo ~ h o 0V^(yf) fe~ Sx2 /° 2 dx 

If Vq(x) = V(x) = \x 4 — \x 2 then yg° = ±1 and yg 1 = =F 4^2 ■ This result shows 
that the positions of the equilibrium states of the central agent will be shifted when 
the individual agents have their own stabilization potential. The states yg - and yg + 
are the two equilibrium states of the central agent, and yg - + (/io/^o)^o( 2 /o _ ) an d 
yg + + (VW(?/o + ) are the two associated equilibrium means of the individual 
agents. 


3. The case of no intrinsic stabilization for the local agents (h = 0) 

In this section we consider the special case where the individual agents have no 
intrinsic stability, i.e., h = 0. In this case, Q is linear so instead of considering 
the empirical distribution T S Xj (t)(dx), we can focus on the empirical mean 

Xjv(t) = jr The pair (xQ N \t) , x n( t)) satisfies the joint SDEs: 

(12) dx^ = —hoV^x^^dt — 0q(xq N) — Xw)dt + ~^= dW°, 

v N 

dxpj ~ —0(x n — Xq N ' l )dt + -^=dWt N \ 

where Wf N> = ^ 7 = ^2f=i W/ is a standard Brownian motion independent of W t °. 

The mean-field limit, (y 0 (t),y(t)) := limAr_ ) . 00 (xg N \t),xw(t)), satisfies |?| with the 
equilibria yg := linp^oo yo(t) = ±1 and y e := linp^oo y(t) = ±1 depending on the 
initial condition (yo(0), y(0)). 
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3.1. Fluctuation analysis in the case h = 0. Here we analyse the fluctuations 
of centred at (yo{t),y(t)) when N is large. To simplify, we assume 

that ?/o(0) = yl = — 1 and y{0 ) = y e = —1, and thus yo(t) = yg = — 1 and 
j/(0) = y e = —1. Define Zq N ^ = VN(xq N ' > — yff) and zn = VN(xn — y e )■ As 
TV —► oo, (zq N \ Zn) converges in distribution to the process (zq, z) where 

(13) dz 0 = -h 0 Vo {yl)z 0 dt - 9 0 {z 0 - z)dt + a 0 dW°, 

dz = — 9(z — Zo)dt + adWt, 


where Wt is a standard Brownian motion independent of Wj) . This means that, 
when N is large, Xg N \t) ~ 2/o + f/N Z o an d x(t) ~ y e + ^7=2 i n distribution. Because 
Vo = y e = ~ 1 ^ the normal state, z o and z are regarded as the central risks (as 
opposed to the large deviations that will be discussed in the next section) of Xg N ' > 
and xj v, respectively. We note that (13) is a system of linear differential equations 
and thus the explicit solution is: 

($)-■“$ 

Therefore (z 0 (t),z(t)) is a Gaussian process with 


)a ( o"o dW°\ 
v crdw s J 


A = 


-hoVoiVo) - do 
6 


(14) 


E (w = " 

(15) 

( Var z 0 {t) 

Co v(z 0 {t),z(t)) 

\Cov(z 0 (t),z(t)) 

Var z(t) 


tA 


zo(0) 

m 


0 (t-s)A ( do 

0 


^- s)A ^ds. 


We want to analyse the impact of the various parameters on (zo(t), z(t)), in 
particular, for the case that t —>• oo and o\ 6 —>• oo with a fixed ratio a := a 2 /6 < oo. 
To do this, we use the eigen-decomposition of A to compute (15) and obtain the 
following. 


Proposition 1. Ifh 0 , 9q and 9 are positive, then\\i\\ t ^ 00 ¥.zg(i) = lim^^ Ez(t) = 
0. In addition, the variances and covariance of the fluctuations zo(t) and z(t) have 
the following limits as t —> oo and a, 9 — > oo with a fixed ratio a = a 2 j9 < oo: 


(16) 


lim 

r ,0 —Voc 

■ 2 /e=o 


lim Var2( 

t—¥ OO 


i (*) = 


2Wfoc ]Y 


(17) 


lim lim Varz(f) 

(7,0 —>00 t—too 

(T 2 /0=Ot 


2 h 0 V''(y e 0 ) 29 ’ 


(18) 


lim lim Cov(zo(t), z(t)) 

(7,6—too t—too 

<7 2 /0—Ot 


2W(% e )' 


This means that after the limits are applied, zq = Z\ and z = Z\ + Zi, where Z\ 
and Z 2 are two independent Gaussian random variables with mean 0 and variances 

a 2 2 

2 fcov 0 ?, (»g ) and res P ectivel y- 

Proof. This involves basic computations given in Appendix | A. 1| □ 
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We see that the variances and the covariance of the limits of zq and z increase 
with increasing cr Q or decreasing h 0 . We also note that these three statistics blow 
up as do — > oo even if <7q/0q is finite and small. This is because when h is exactly 
zero, Xn cannot serve as a stabilizing term and Xg N ' > cannot diversify its risk to Xn 
by increasing 6g. 

3.2. Large deviations. 

3.2.1. A general large deviation principle. From the mean field and fluctuation 
analysis we see that if N is large and x^\o) = 27,(0) = —1 for all j = 1 ,... ,N, 
then one can expect that (xg N \t), x n( t)) ~ (yg,y e ) = (—1, —1) for all t. However, 
as long as N is finite, Xg\t) and xjv(f) are stochastic processes and therefore 
the event that the overall system has a transition in a finite time interval has a 
small but nonzero probability. Mathematically speaking, we consider the event of 
the continuous paths (xg N \t), x n (t)) £ C([0, T], R 2 ) starting from {yg~,y e ~) := 
(—1, —1) at time 0 to ending around (jjg + ,y e+ ) := (1,1) at time T: 

(19) A = {(a;o(t),2’(f))te[o,T] € C([0,T],R 2 ) : 

(so(0),®(0)) = (-l,-l),||(*o(T),*(r))-(1,1)|| <<5}, 

where || • || is the standard Euclidean norm in R 2 . 

The Freidlin-Wentzell theory [8) Section 5.6] says that, for N large, P((xq N \xn) £ 
As) satisfies the following large deviation principle: 

- inf I{x)< lim inf logP((j;[ ) Ar) , x N ) £ As) 
xeAs n^oo 

< lim sup logP((a;^ Ar) ,5iv) € As) < - inf I(x), 

IV—>oo -'V x£As 

where As and As are the interior and closure of As under the standard C([0, T ], R 2 )- 
topology, respectively, and /( x) is the rate function for the exponential decay of 
the probability that will be specified later. By using a similar argument as in PH 
Lemma 5.2], we can show that for any e > 0, there exists sufficiently small S > 0 
such that 

— inf I{x) < liminf logP((2’i Ar \ Xn) £ As) 

x^A N—yoo 

< limsup -j- logP((xg Ar \ xn) £ M5) < - inf J( x) + e, 

A—>00 Al x£A 

where 

(20) A= {(x 0 (t),x(t)) te[0tT] £ C([0,T],R 2 ) : 

(®o(0), *(0)) = (-1,-1), (® 0 (T),x(T)) = (1,1)}. 

In other words, for large N and small S , 


and we define this probability as the systemic risk of the overall system. We will 
discuss the rate function I(x) separately for the cases that a 0 = 0 and cr 0 > 0 in 


( 21 ) 


P(( a '0 JV ^, x n) e As) « exp 


—N inf I(x) 
xG A 
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the following sections. We will next compute the minimum of the rate function 
infxe.4/(x) to obtain the systemic risk in (21). 


The minimizer x* = arg min ;/ .p_4. l{x) is the most probable path for the 
rare event As in the sense that the mass of the conditional probability P(-|^) is 
concentrated around x* exponentially fast as N -A oo . Indeed, if x* exists and is 
unique, then for any open neighbourhood N(x*) containing x*, 

(22) ¥((x^ 0 N \x N ) e N(x*)| (x { 0 N \x N ) £ As) 

= 1 - P((xq N \xn) £ 1SS (x*)\(xq N \ x n) £ As) 

P((x { 0 n \x n ) eN c (x*)n As) 


= 1 - 


P((x[, JV) ,a;jv) £ As) 

> _ exp(-iVinf a; g N c (a; , )n ^ i5 1(x)) N—>oo ^ 

exp(-NM xeAs I(x)) 


by using the fact that x* is unique and As is closed. 


3.2.2. Degenerate case. We first consider the degenerate case where do = 0 and 
er > 0. Then (12) becomes 


JLW - -h n V'(x m 
^ X 0 — ^O^OV^O 


) - 9o(x ( 0 N) - x N ), 


dxN = —9{xn — Xq N ' > )dtt H— -z=dWf N K 

v N 


The rate function I(x) in (21) is of the form 

1 f T 

(23) I(x) = I(x 0 ,x) = ^2 J (x(t) + 9(x(t) - x 0 (t)) 2 dt, 

if (x(t))t£[o,T] is absolutely continuous in time and xq = —hoVg(xo) — 9q(xq — x) 
and I(xq,x) = +oo otherwise. Here the dot stands for a time derivative. By (21), 
in order to compute the systemic risk, we need to solve the optimization problem: 


(24) 


1 f T 

inf tt-o / (S(t) + 9(x(t) - x 0 (t)) dt, 

x(t) 2 <t z J q 


with the constraints that (£(t))t<=[o,T] is absolutely continuous in time, xq = —hoVg(xo) 
9o(xo ~ x), Xg(0 ) = 5(0) = —1 and xq(T) = x(T) = 1. By using x = ^d’o + 
^V'(xg) + Xq, the constrained optimization problem is equivalent to 


(25) inf 

Xq Z(J 


7TX o + ^-Vq(x 0 )x 0 + (1 + ^-)xo + 6 -^-Vq(x 0 ) 
Vq V 0 v 0 u 0 


dt, 


with the boundary conditions xo(0) = —1, Xq(T) = 1 and i’o(0) = Xq(T) = 0. From 
basic calculus of variations, the minimizer Xq satisfies a fourth-order boundary value 
problem that we describe in the Slowing proposition. 
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Proposition 2. The minimizer (xq,x) of inf ( X0jS ) £- 4 I{xg, x) of the rate function 
p?3| ) satisfies the following boundary value problem 

(26) 


d 4 . 2 d 2 

^x 0 -(9 0 + 6) —x 0 


VI'" (x 0 ) 


- 9oVo"(x 0 ) ( -jjX 0 ) - 29qVq'{xq) 


dt 

d?_ 

dt * 2 


x 0 ) + 3Vq'(x 0 ) 


dt 


x 0 


dt 2 


x 0 


x 0 


+ h 2 0 V'\x 0 ) 


~K"( x o) {jff t x ^j - n'(x 0 ) (Jj^x + 0 2 Vq{x q ) 


= 0 , 


with Zo(0) = -1, x 0 (T) = 1, £x 0 (0) = fjfx 0 (T) = 0, and 

x(t) = + ^V'{x 0 (t))+x 0 (t). 

Proof. See Appendix |A.2| 


□ 


If ho = 0, we can solve xo and x explicitly. The boundary value problem (26) is 
then 


(27) 




with the boundary conditions x o (0) = —1, fjfx 0 ( 0) = 0, x 0 (T) = 1 and J \xo{T ) = 0. 
The associated minimizer x is x(t) = x 0 (t) + f-4iXo{t). The solution of (27) is 


(28) x 0 (t) = 

(29) x(t) = x 0 (t) 


(1 + e- (eo+e)T )(2t-T) + 


do dt." 

2 -(G n +6)t __ 2 (8 n+ e)(T-t) 

(e 0 +e) e (e 0 +e) 


T( 1 + e-( 0 o + e) T) + _JL_ {e - { e 0+S)T _ i) 

2 (1 + e ~(®°+®) r ) — g — (®o+0)t _ g— (#o+0)(T—*) 


0 O T(1 + e-( 0 o+«) T ) + jg^e)(e-^+ e ) T ~ l) ' 

These are the most probable paths followed by the two processes to realize the 
rare event asociated with the systemic risk. Note that x(t) is ahead of xg (t), which 
means that the individual agents drive the transition. We also obtain the following 
proposition. 


Proposition 3. If ho = h = 0, then the probability of transition is 
(30) 


P((x ( 0 N \x N ) e As) 


IN(9 0 + Of 




For large T (i.e. (0 q + 0)T 1), the most probable paths are 

2 1 

(31) x 0 (t) ~ x(t) « -1 + —, 
and the probability of transition is 

(32) P((x { 0 N \x n ) g As) » exp ) ■ 

This shows that stability increases with 6 and decreases with 0 O - This is because 
when <To = 0 and cr > 0, xg is a stabilizing term while x is a destabilizing term. 
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When 9 increases, x (unstable) is forced to be close to xo (stable), and therefore 
the systemic risk is reduced. On the other hand, the systemic risk is higher if 9q 
increases, as we make Xq stay close to x. 


3.2.3. Non-degenerate case. We next consider the non-degenerate case where <jq 
and cr are positive. In this case, the rate function I(x) in (211 has the form 
(33) 

1 r T 1 f T . 

I(x) = I{x 0 ,x) = ^2 J (xo+h 0 Vo(x 0 )+9 0 (x 0 -x))' 2 dt+ 7 ^ J (fi+0(x-x o )) 2 dt, 


if (xo(t))tg[o,T] and (S(t))tg[o,T] are absolutely continuous in time and I(xo,x) = 

+00 otherwise. Again by the calculus of variations, the minimizer (x 0 , x) of inf( Xo 2 ) e ^ J(xo, x) 
satisfies a system of second-order ordinary differential equations. 


Proposition 4. The minimizer ( Xq,x ) of inf( X0)S ) £- 4 I(xq, x) of the rate function 
(33) satisfies the following system of second order boundary value problems 


( 34 ) ^ x 0 = ~^(<? 2 9 o - a 2 9)j t x + ^( ct 2 6 ^ + &o0 2 )(xo - x) 

+ h 0 9 0 [Vo(xo) + Vo'{x 0 ){x 0 - x)] + hoV r 0 , (x 0 )V r 0 "(x 0 ) 

= ~ 2 ^ a o 9 - v 2 8o)^_x 0 + ^(o'o 6 ' 2 + o" 2 6»o)(x - x 0 ) - h 0 




with Xo(0) = x(0) = —1 and Xq (T) = x(T) = 1. 

Proof. The proof is essentially the same as the proof of Proposition [2] in Appendix 
IA.2I and thus is omitted. □ 


Although (34) is solvable when ho = 0, the explicit solution is very complicated 
even for zero ho- Therefore we compute the transition probability by using the fact 
that (xo (T), x(T)) are jointly Gaussian random variables and obtain the exponential 
rate of the decay of the probability. 


Proposition 5. If ho = h = 0 and Xo(0) = x(0) = — 1, then the probability of 
transition has the following exponential rate of decay: 


(35) 


P((xo Ar) ,Xiv) e As) « exp N 


2(9 0 + 9) 2 


T(9 2 Oq 


9 2 o<x 2 ) 


for large T. 

Proof. See Appendix |A.3| 


□ 


3.2.4. The case that ho > 0. Most of the large deviation analysis in this section is 
about the case ho = 0 in order to have explicit results. Although it is also possible 
to consider the case that 0 < ho 1 and use the small ho analysis, we will solve the 
large deviation problems numerically as the associated boundary value problems ([2]) 
and (34) can be solved easily by standard numerical methods. The details of the 


numerical analysis are presented in Section [6j 
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4. Formal large deviations for the empirical measures 


In this section, we extend the large deviations formulation from the space of real¬ 
valued processes (xQ N \t), £iv(i))te[o,T] to the space of probability-measure-valued 

processes (xo N \t),U N (t,dx)) te [ 0tT ], where U N (t,dx ) := ^ SyLi d Xj ( t ){dx). The 
reason we consider a more general and complicated space is that there is no closed 
equation for xn when h > 0, because Q is not linear for non-zero h. In addition, we 
obtain more information by considering the more general space even for h = 0 and 
we show that when h = 0 the generalized problem is (at least formally) equivalent 
to the problem we considered in the previous section. 

We also note that there are no existing large deviation results for (x^\t), C/jv( t, da : ))te[o,T] 
satisfying ([3]) and ([4]) even if h = 0; the current most general large deviation princi¬ 
ple for weakly interacting particle systems is [S], but unfortunately our model still 
cannot be covered. Thus the results in this section are formal. 

Motivated by [ 6 ], the (formal) rate function for (xQ N \t), f7jv(t, dx)) te [o,T] satis¬ 
fying © and Q is 

1 fT 

J{(xo(t),(j>(t,dx)) te[0t T]) = / (xo + h 0 Vo(x 0 ) + e 0 (x 0 -x)) 2 dt 

Z(J 0 Jo 

1 f T (<t>t - h-^[V'(x)4>] - \v 2 <t> xx - 0^[(x - x o (t))0],/(a ;)) 2 

+ ^L - mnm - dt • 

for er 0 > 0 and for ctq = 0 , 


J({x 0 (t),<t>(t,dx)) 

te[o,T]) 


1 

2 ct 2 


sup 

f(.x):{c/>,{f'(x)) 2 )^ o 


- h -LW o*#] - \° 2( t>xx - d-jL[{.x - x 0 (t))^],f(x)) 




(37) 


p{t,x) = 


1 


2 ^ 


: exp 


(x - x{t)Y 

2 — 
z 28 


-dt , 


if xq + hoV^xo) + do{xo — x) = 0 or J({x 0 (t), <j>(t, cte))te[o ,T]) = oo otherwise. Here 
/ is in the Schwartz space, (</>,/( x)} = / f(x)<p(t,dx), and the partial derivatives 
(Jj, ^ 2 ) are defined in the weak sense. 

By the contraction principle [5] Theorem 4.2.1], if the large deviation princi¬ 
ple for (xq N \t),Uisr(t, dx)) tB [ 0 T] exists, then by using the projection XQ N \t) <—> 
x o N \t) and Uiy(t,dx) 1 —> Xjv(t) = {UN(t,dx),x), the large deviation principle for 
(xq N ^ (t) , xjv(t))te[o t] a ls° exists with rate function 
(36) 

The following result shows that when h = 0, for either er 0 = 0 or a 0 > 0, 
Z((xo(t),x(t)) te[0 ,T]) = I , x(t)) te[0>T] ) in or d33h, respectively. 


Proposition 6 . If h = 0, then the infimum in (36) is reached for and only for the 
path of Gaussian density functions 


In addition , l((x 0 {t), x(t)) t& [ 0 ,T}) = l({x 0 (t),x(t))t^T]) in (23) for a 0 = 0 and 
in p5ll for cr 0 > 0 . 
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Proof. See Appendix [B| 


□ 


In other words, when h = 0, we can simply consider the large deviation problem 
for (xo(i), xjv(i)) in Section [3] instead of (xo (t),Xjsr(t,dx)) in a complicated space. 

However, if h > 0, then it is necessary to consider (x^ (t), U]y(t , dx))t^[o,T] with 

rate function J ((xo(t), 4>(t, dx)) t ^[o,T]) as now the large deviations for (xQ N \t),XN(t, dx)) te ro j T] 
cannot be obtained by the Freidlin-Wentzell theory. Motivated from Proposition 
[ 6 ] and QB Section 7], we know that because for h = 0, the most probable path 
for the empirical measure UN{t,dx ) is the Gaussian probability measure p(t,x)dx, 
it is reasonable to assume that for 0 < h -C 1 , the most probable [/jv(t, dx) is a 
Gaussian probability measure plus higher order corrections in h. In addition, as 
the base case [h = 0) is Gaussian, we parametrize the most probable path of the 
density </>(i, x) by the Hermite expansion: <f> = p + hq W + /i 2 qG) + • • •, where 


P(t,x) = 2. - exp , p{t) = {(j){t,x)dx,x), 


2 ’rfe 


J 2 e 


oo 077, 077, 

q (1) {t,x) = ^2/3 n (t)—p(t,x), q (2) {t,x) = '^2'y n {t)—p(t,x). 


n—2 


n—2 


Then 


(j)(t, dx)) te[0 , T ]) = min J ((x 0 (t), /x(t), Pn(t), 7n(0)te[o,T])+o(/i 2 ), 
Xo,/I,Pn,7n 

and we can solve the associated variational problems for Xq (£), fi(t), /3 n (t) and 7 n (t) 
as in pn Section 7]. This task is not carried out in this paper. 


5. Optimal control of the central agent 

In this section, we consider an optimal control problem by introducing a control 
term a.j{i) into Q. In order to be able to address the problem in a manageable 
way and to discuss the role of the parameters, we will write it as a linear-quadratic- 
Gaussian control problem as in [3]. We let h = 0 and define X^ ( t ) = x^ (£)—j/g = 
Xq N) (t) + 1 and Xj(t) = Xj(t) — y e = Xj(t) + 1. By assuming that Xg N \t) is small 
so that hoVo(xQ N \t)) = hoV^y^ + X^ N \t)) « H 0 X^ N \t) with H 0 > 0, we have 

1 N 

(38) dX ( 0 N) = -H 0 X { 0 N) dt - 9 0 (X { 0 n) - X N )dt + X N = ~Y, X i> 

V A iV j= i 

(39) dXj = —0(Xj - X ( 0 N] )dt + crdWf + oijdt, j = 1,..., N. 

The optimal controls ay are adapted to the past {(X,(s))j = o,...,jv>0 < s < t} and 
such that the following cost function is minimized: 

1 N 

(40) J(ai,...,aiv) = 

3 = 1 

This cost function means that the optimal controls try to make X,- close to X^ N> 
with a quadratic cost. We can regard the term —0(Xj—X^ N ' > ) as a passive feedback 
while ctj is the active feedback from the central agent. A possible control (but not 
optimal as we will see) is to take the active feedback ctj = —9 c (Xj — Xq N ^) for 


a 2 At) + e 2 c {X { 0 N \t)-X 0 {t)) 2 dt 
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some well chosen 0 C . The goal of this section is to study the form of feedback that 
the optimal control produces and whether it is different from the passive feedback 
— 0(Xj — X ( j ,V ^). By using standard theory, we have the following optimal control 
aj{t) for (X ( 0 N \t),X N (t)). 


Proposition 7. The optimal control etj(t) that minimizes J in (40) where (Xq N \ t), Xjv(t))tg[o,T] 
satisfies \3tfy and (3fJ\ ) is 

(41) aijft) = -0 C (b(t)X^ N) (t) +d(t)Xj(t) + e(t)X N (t)^j , j = 

where (a(t),b(t),d(t),e(t)) tG [ 0 ,T] tde solution of the following Riccati equations: 

(42) 

a(t) = 2 (0 O + H 0 )a(t ) — 2 06(f) + 9 c b 2 (t) — Q c , 

b(t) = ( 0 o + H 0 + 9)b(t ) — 9d(t) — 9oa(t) + 9 c b(t)d(t) + 0 C — 0e(t) + 9 c b(t)e(t), 
d(t) = 2 9d(t) + 0 c d 2 (t) — 0 C , 
e(t) = —2 9 0 b(t) + 2 9 eft) + 9 c (2d(t)e(t) + e 2 (t)), 
with the terminal conditions (a(T),b(T),d(T),e(T)) = (0, 0,0, 0). 

Proof. See Appendix [Cj □ 

When T —>■ oo we have 


(43) 


-j(t) = ~°C ( 


boo X () 


(A) 


( t ) + doaXjft) + eooX^v 


■ / C L, OO u 'CX)J 


where the parameters (cioo, &oo, doo, e^) satisfy the algebraic Riccati equations: 
(44) 0 = 2(9 0 + H 0 )a oo -29b oo +9 c bl o -0 c , 

0 = (0o + H 0 + 9)b oo — 0^ — 0o a oo + Ocboodoo + 9 C — Oe + 9 c b 0 
0 = Odoo + 9 c d 2 0 — 0 C , 

0 = —20q0oo + 20COO + 0 c (2d oo e oo + e^). 

In these conditions (Xq N \Xn) satisfies the SDE: 


dX, 


(AT) 


0 = -HoX^’dt - 0 o (X^> - X N )dt + 


dXjy = -0(Abv - + -JLdW t (JV) - 9 C 

v iV 


where Sj=i W^t) is a standard Brownian motion. 

In order to obtain the optimal control (43), we need to have the coefficients 
(boo, doo, eoo) that cannot be obtained analytically, in general, and must be com¬ 
puted numerically. However, we are able to find approximate solutions in certain 
regimes. We note that from (44), doo = (—9 + \JO’ 1 + 9"*)/9 c , and we consider the 
following cases: 

(1) If 0o = 0 and Hq = 0, then we find boo = —doo and eoo = 0, so that we 
obtain the system 


L o 


boo X o 


(AT) 


(do 


dt , 


dX { o N) = -HoX^’dt 


(N), 


^dWl 

Vn * 


dX N = -?=dW t {N) - y/0 2 +9 2 (X N - Xo N) )dt, 

v N 
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which shows that the passive control —9(Xj—X q N ^) and the optimal control 
ctj combine in a quadratic way to form the feedback — \/9 2 + 9 2 (Xn — 

4 N) )- r _ 

(2) If 0 < 9 0 <C 1 and H 0 = 0, then we find b^ = — d oo +9 0 d oo /\/9 2 + 9 2 +o(9 0 ) 
and eoo = -9 0 d tx /^/9 2 + 91 + o(9 0 ), so that we obtain the system 


dX. 


(N) 


o = -0„( X { 0 N) - X N )dt + -^dW t u , 


CTO 


- X { 0 N) )dt 


dX »- -h dWim - (^ - 9 “^Sr) 

which shows that the optimal control chooses to reduce the feedback, prob¬ 
ably because is destabilized by 9 0 . 

(3) If 0 < 0o ^ 1 and 0 < Hq <C 1, then we find = —doo + {Hq + 
do)doo/x/0 2 + 91 + o{9 0 ,H 0 ) and e x = -9 0 d oo /^/9 2 + 9 2 + o(9 0 ,H 0 ), so 
that we obtain the system 


dX ( 0 N) = -H 0 X^’dt - 9 0 (X^’ - X N )dt + 


AN) 


(N) 


oo 


>-0 — ^ is)™ i 


dX N = -°=dW t (N) 

Vn 


-Hr 


y/9 2 + 9 2 -9 


-X ( 0 N) )dt 


-Xjsrdt, 


y/9 2 + 9 2 

which shows that the optimal control chooses to reduce the feedback but it 
also controls Xn directly. 

6. Numerical results 

6.1. Numerical results of fluctuations. In this subsection we compare the ana¬ 
lytical fluctuation results ( I6|l8 ) with the fluctuations obtained from the numerical 
simulations of (xq N ^ (t), Xw(t)) in (12). We use the Euler scheme to discretize (12): 


(45) 


Xo N) (n+ 1) = -^=AW° +1 - h 0 Vo(x { 0 N \n))At - 9 0 (x { 0 N \n) - x N (n))At, 


x N (n + 1) = -j=AW n+1 - 0(x N (n) - x, 
AN) 


(A) 

0 


(n))At, 


with Xq (0) = a;jv(0) = —1 and {A W° +1 } n , {AW n +i}„ i.i.d. Gaussian random 


variables with mean 0 and variance At. We simulate (45) up to time T and we 
take T large enough so that (xQ N \t),XN(t)) is in equilibrium after T/10. Therefore, 
Var(lim t _ ) . 00 Var(lim t ^ 00 x N (t)) and Cov(lim t _ ) . 00 x^ N \i), lim^oo x N (t)) 

are approximately the sample variances and sample covariance of {xq N ^ ( n ) : T/10 < 
nAt < T} and {x N (n) : T/10 < nAt < T}, respectively. 

For each simulation, we vary one parameter for 100 different values equally dis¬ 
tributed in the region of interest, and use the values in Table [Tlfor the other pa¬ 
rameters. The results are shown in Figures |T] and [2j In Figure]!] we compare the 
analytical formulas ( 16||18| ) with the sample variances and sample covariances from 
the direct numerical simulations for 100 different ho and <To uniformly distributed in 
the region of interest. In Figure [2] we compare the analytical formulas (T6|[T8) with 
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N 

T 

At 

ho 


0o 

a 

9 

100 

10 3 

10" 3 

0.5 

0.1 

0.1 

1.0 

10 


Table 1. The typical values of parameters used in Sec |6.1| For 
each simulation, we vary one parameter and the other parameters 
are fixed at the values in the table. 


3 

2 

0 

3 

2 

0 

3 

2 

0 



Variance of Fluctuations of Xq x 10 -3 Variance of Fluctuations of Xo 


- numerical 

4 

- numerical 

— — - analytical 


| - - - analytical | 


2 


0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 0 

Variance of Fluctuations of xn 

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 

10 -3 Variance of Fluctuations of xn 




0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 

j-3 Covariance of Fluctuations of Xq and x n 



0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 

h 0 



0 1 1 '- 1 - 1 - 1 - 1 - 1 - 1 - 1 - 

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 

x iq- 3 Covariance of Fluctuations of Xq and x n 



Figure 1. We compare the analytical formulas for variances and 
covariances with direct numerical simulations. On the left the hor¬ 
izontal axis is ho and on the right a q. 



e 


Figure 2. Same and in Figure [l] except that the horizontal axis 
on the left is a and on the right 9. 


the sample variances and sample covariances from the direct numerical simulations 
for 100 different a and 9 uniformly distributed in the region of interest. We see that 
there is goo d agreement between the analytical formulas and the simulations and 
thus (16 


18) indeed capture the fluctuations of the equilibrium of (xg N ^ (t), xjv(t)). 


6.2. Numerical results of large deviations. In this subsection, we compute 


the most probable paths ( Xo,x ), defined in Section 3.2 by numerically solving 
the associated boundary value problems (26) and (34) for cr 0 = 0 and a 0 > 0, 
respectively. We use the boundary value problem solver bvp4c in MATLAB to 
solve these problems. The details of the algorithm can be found in [T5j . 
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For the non-singular cases, for ho small, we use xq( t) = — 1 or xq (t) = (2 t/T) — 1 
for (26), and x 0 (t) = x(t) = —1 or x 0 (t) = x(t) = (2t/T) — 1 for ([34]), depending 


on which one gives better results. We found that bvp4c sometimes did not give an 
accurate solution even for the non-singular cases. The numerical solutions failed 
to pass their internal accuracy check of the MATLAB routine. The reason for this 
is not clear. However, this issue can be bypassed by iterating bvp4c several times. 
More precisely, we use the inaccurate solution as a new initial guess and use bvp4c 
to solve the same boundary value problem again to obtain a new solution and so on. 
After several iterations, bvp4c finds the correct solution that passes its accuracy 
check. 

For the nearly-singular case, when ho is large, the method just described fails 
to find the correct solutions even with several iterations. To get past this issue, 
we use as initial guesses solutions of the less singular cases obtained by the above 
technique. For example, we use the solution of the problem with ho = 1 as an 
initial guess to solve the problem with ho = 2, and so on. Eventually we can solve 
some quite singular problems, for example, with ho = 10 . 


6.2.1. Impact of ho . In Figure[3]we plot the most probable paths (ly, x) as functions 
of time, for ho from 0 to 10. On the left all the plots are with ctq = 0 and on the 
right (T 0 = 0.5. We note that when h 0 = 0, (xq,x) is smooth and in fact it is 
approximately linear, while (xo,x) is quite curved for h 0 = 10. We see that when 
Xo(t) < 0, the destabilization of the system is driven by x(t). Indeed, x has higher 
external risk (a = 1) than xq( t) does (< 7 o = 0 or ao = 0.5) and has no intrinsic 
stability (h = 0), and therefore in the most probable path x(t) destabilizes xq( t). 
Nevertheless, once Xq( t) > 0, the system transition is driven by Xq( t) because the 
double-well potential forces x 0 to go to the failed state 1 , and x(t) is driven by x 0 (t). 
This effect is strengthened when ho is large because the double-well potential plays 
a more important role in that case. 

In Figure [4] we plot the values of inf x 6 ^/(a;) for different ho . We see that 
inf a. e _ 4 1(x) is an increasing function of ho- This is expected because the system 
is more stable if it has more intrinsic stability (ho)- We also see in Figure [4] that 
bif^e-A I( x ) has quadratic behavior with respect to ho for small ho and linear be¬ 
havior for large ho- 


6.2.2. Comparison between small fluctuations and large deviations. Here we com¬ 
pare the small fluctuations of (xq N ^, xn) described by the processes Zo and z in (13) 
and the large deviations of (xq,Xn) described by the infimum of the rate func¬ 
tion mi x ^I(x). For the characterization of the small fluctuations, we compute 
limt-xx, Varzo(i) in (49) and lim^oo Varz(f) in (50). For the characterization of 
the large deviations, we compute I(xo,x) in (23) for ao = 0 where (xq,x) is the 
solution of (26) and compute I(xo,x) in (33) for ao = 0.5 where (a;o,x) is the solu¬ 
tion of (34). The goal is to visualize the fact that the systemic risk characterized 
by infa , g _4 I(x) may vary significantly even though the individual risk measured by 
lim^oo Varz(f) is kept at a fixed level. 

Motivated by ( [16]) and ( p~7[ ) , we know that linp^oo Varz 0 (t) and lim^oo Varz(t) 
are not significantly affected if we increase a and 9 but keep the ratio a 2 /9 the 
same. In Figure[5]we confirm this expectation and we also observe that inf^g^ I(x) 
increases as a increases, which means that systemic risk decreases. This also means 
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I = 0.88899, T = 10, Ao = O,0 o = 1,8 = 1 ,<t 0 = 0,<r = 1 



t 

I = 1.2943, T= 10, h 0 = 1,9 0 = 1,9 = l,<r 0 = 0,ir = 1 



t 



t 



t 



Figure 3. The most probable paths (xo, x) = arg min _4 I for h o = 
0,1, 5,10. We let T = 10, 6 0 = 1, 9 = 1 and a = 1. The left column 
is the case cto = 0 and the right column is the case er 0 — 0.5. 
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Figure 4. The infimum of I over A: inf _4 / for ho = 
0, 0.1,0.2,..., 1 and for h 0 = 0, 1,2,..., 10. We let T = 10, 9 0 = 1, 
9 = 1 and a = 1. The left column is the case cr 0 = 0 and the right 
column is the case <7o = 0.5. 


that, for a fixed level a 2 /9 of individual risk, the reduction of 9, ie the interaction 
of the local agent with the central agent, reduces the systemic risk. 

One may also expect that 9 0 does not greatly affect lim t _,.oo Varz 0 (f) and lim^oo Var5(f) 
however, in Figure[ 6 ]we see that the effect of 9q on lim^oo Varz 0 (f) and lim^oo Var z(t) 
is not negligible. In other words, the independence of lim t _>.oo Varz 0 (f) and lim^oo Var z(t) 
with respect to 9 0 only holds in the limits (16) and (17). 


6.3. Numerical results for optimal controls. In this subsection, we use the 
Euler scheme to simulate © with optimal controls: 

(46) XQ N \n + 1 ) = -^AW° +1 -hoVo(x i 0 N \n))At-9o(x i 0 N \n)-x N (n))At, 

x N (n + 1) = — —AW n+ i — 9{xn{ti) — x\^\n))At + a°°{n)At 
V N 

with £ 0 ^( 0 ) = 5jv(0) = —1 and {AW° +1 } n: {AW n+ i} n i.i.d. Gaussian random 
variables with mean 0 and variance At, where 

(47) af (t) = -9 C (b oo {x ( 0 N \n) + 1) + d 00 (x j (n) + 1) + eoo^ivW + 1)) 


and ( 1100 , 600 ,^ 00 , 600 ) satisfies the algebraic Riccati equations (44). 
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Figure 5 . Plots of mi x ^I{x), lim^oo Va.rzo(t) and 
linit->oo Varz(f) for cr 2 from 1 to 10 with cr 2 /9 = 1 . We let 
T = 10 , 9q = 1 . The left column is the case do = 0 and the right 
column is the case <j 0 = 0 . 5 . 


To obtain (aoo, boo, doo , e-oo ), we numerically solve (42) for large enough T so that 
(a(0), 6(0), d(0), e(0)) is essentially (aoo,boo,doo, eoo). The values of the parameters 
used in (46) are listed in Table [2j 

We see from Figure [7] that the uncontrolled problem is very unstable in the sense 
that XqN' 1 and Xn jump frequently between —1 and +1. On the other hand, under 
the same values of the parameters, the controlled and are much more stable 
with no transition from —1 to +1. 
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T= 10,ft o = 1,8 = 10, ito = 0,<t = 1 T = 10, /ip = 1,0 = 10, <r 0 = 0.5, a = 1 






Figure 6. Plots of mi x ^I{x), lim^oo Varzo(i) and 
lim^oo Var z(t) for 9q from 1 to 50. We let T = 10, 9 = 10, a = 1. 
The left column is the case <7o = 0 and the right column is the 
case cr 0 = 0.5. 
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Figure 7. Sample paths of x^ N \t) and a;jv(t) with and without 
the optimal control. With the optimal control, XQ N \t) and xi v(t) 
are much more stable than the uncontrolled ones. 

7. Summary and Conclusions 

We have formulated and analyzed a multi-agent model for the evolution of indi¬ 
vidual and systemic risk when there is a central agent acting as a stabilizer in the 
system. The local agents do not have an intrisinc stabilizing mechanism. The main 
result of this paper can be visualized in Figures [5] and [6] and is briefly described 
as follows. The systemic risk decreases when the rate of adherence of the local 
agents to the central agent increases, but it increases when the rate of adherence 
of the central agent to the mean of the local agents increases. This is under the 
condition that the observed individual risk is kept approximately constant. We also 
show that the effect of drift controls on the local agents is to always stabilize the 
systemic risk. 
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Appendix A. Proofs in Section [3] 


A.l. Proof of Proposition [l] We first consider the eigen-decomposition of A: 
A = QAQ -1 , where 


A = 


Ai 0 

0 A 2 


Q = 


Ai — A 2 


i + 4 i + 4 


Q = 


i -(i + 4) 


-i 


1 4- 
1 ^ e 


*1 = 1 {-[W4o e ) + 00 + 0} + 0/ioW) + 00 + or- - 40fco^'fe§)} , 
A 2 = \ {4W(z/o e ) + do + 0} - sJ[h 0 V£(y e 0 ) + 9 0 + 0] 2 - 49h 0 V''(y e 0 )\. 


We note that Ai and A 2 are real and negative if h 0 , 9q and 9 are positive. Then from 
(fl4|, lim^oo Zo(t) = limj^oo z(t) = 0. In addition, from the eigen-decomposition 
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we have 
(48) 


f Varzo(i) Cov(z 0 (i), z(t)) 
yCov(zo(t),z(t)) Var z{t) 

t 

= Q ' 


j o e (t ~ s)A Q 1 (o“ ° 2 ) (Q- 1 ) T e (t - s)A rf S Q T . 


We observe that 

t 2 

0 a 2 


Q 


-1 / °0 0 


(Q 


- 1 \T 


-o 2 + - 2 (l +¥) 2 


_<J o - cr2 (l + ^)(1 + if) °o + <j2 (l + T 1 ) 


-a 2 -a 2 (l + ^)(l+^) 


0 L/ V x T- -g-;v x 1 IT 

g + a 2 (l + ^ 2 


Then 


aTTWI - _°o - ^ 2 (1 + if )(1 + W 

0 


-2X7^0 + tj2 ( 1 + if) 2 ] 

1 r ^.2 _ _l Ai 
e 


A1+A2 [~°0 “ Cr2 ^ 1 + if H 1 + if)] ~2Xl[°0 + <j2 ( 1 + if) 2 ] 


So we obtain 


(49) 


t limVarz 0 (t)= (X —^ — 


1 + 


Ai + A 2 

1 

~2\ 2 


1 + ^ I ( 1 + ^ 


(Jf\ (J [ 1 + 


Ai 


(j n +cr 1 + 1 + - 7 - 


1 + v 


(50) 


(Jr, + (7 ( 1 H-— 


1 


hm Varz(t) = ^ x | erg + a 2 (1 + ^ 


Ai + A 2 


(Ai - A 2 ) 2 I 2Xi 


cr n+ cr “(l + — ) (1 + wr 


1 

2 W 


<Tn + O 1 H- ~ 


( 51 ) 


t hm Cov(z 0 (t),z(t)) = l 1 + f 


+ 

+ 


1 


Ai + A 2 

1 

Ai + A 2 

1 


1 + 


e 


1 + v 


Ai 

e 

A 2 

T 

Un + cr“ ( 1 + — ) ( 1 + — 


a 0 2 +a 2 (l + ^ 


A 2 

~e 


<Xg + <T ( 1 + 


Ai 


1 + 


2A 2 


1 + 


A 5 


CTq + (7 ( 1 + 


Ai' 
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We are interested in the case that er and 9 go to infinity while the ratio a = a 2 /9 
is fixed. For 9 large and using the approximation y/1 + x = 1 + \x + 0(x 2 ), we 
have the following expansions: 


A] 


29 


-[«(% e ) + + 0] + [hoWo) + + 9} \ /l — 


hoK(y e 0 ) 


4 9h 0 V’’{y^) 


{h 0 V£(y e 0 ) + 9 0 + 9] 2 


0,1 

h 0 V''(y e 0 )+9 0 +9 ' W 2 


A' 


1 


1 + f = 29 ' W ~ [WfoS> + ^0 + 0] - [^o^o"(2/o) + ^0 + 0] A /1 - 


4^oF 0 "(yg) 


— ~a[hoVo(yo) + Qq] + 


lh 0 V£(y e 0 ) + 9 0 + 9} 2 


hoVoiVo) 


O 


9 2 


hoVoiUo) + 9 q + 9 

Thus Ai —¥ h 0 Vo(yg) as 9 —> 00 and 1 + ^ = O(^) and finally we have the limits 
and 


A.2. Proof of Proposition [2} If Xo is the minimizer, then for any perturbation (f> 
with ^>(0) = 4>{T) = </>(0) = </>(T) = 0, the directional derivative of I must be zero: 


d 

de 


1 f 1 

H*0+e<l>) = ^J 


+ + (f + io + 


e—0 

7 T 0 + 7 T V o'( x o)<l>xo + ] j-Vq (x 0 )j> + (l + 1") </> + ^-Vq(xo)4> 
(to (to (to v ‘A) / FQ 


'0 

dt = 0. 


After integration by parts and using the fact that (f> is arbitrary, the minimizer Xq 
must satisfy the following equation: 

9 0 dt 2 

l^o + ^-V^(x 0 )x 0 + ( 1 + 1) xo + ^-Vo{x 0 ) 

Vo ^0 V Vo J Vo 


J-Xo + Y V o(xo)x 0 + (1 + 1) X 0 + ^-Vo(xo) 
Vo v 0 \ v 0 / v 0 


+7tK"( x o) x o 

(to 


d ( h 0 


dt ^ 0o 

-fi + 


VS{x 0 ) 


+ TT^tAcO^o + f 1 + io + -^v'(x 0 ) 

(to (to \ (to / (to 

1 h 0 a, . ( 9 \ 9ho i, . 

x-x 0 + v-P 0 (x 0 )a;o + 1 + 7“ P’o + -^-^(* 0 ) 
(to (to V (to / (to 

7^0 + -n-Vo{xa)xa + (l + — 'j x 0 + ~w—Vo( x o) 

(to (to V (to / (to 


d 

7q J dt 


1 T/■///•_ \ 

^0 ( x o) 
(to 


= 0 . 


with the boundary conditions xo ( 0 ) = — 1 , xo(i) = 1 and -^xq ( 0 ) = -^xq( t) = 0 . 
We then obtain (f26j) after rearranging the above equation. 


A.3. Proof of Proposition [5j If /iq = 0, (12) is a system of linear SDEs, and the 
explicit solution can be found: 


ra - 

\x n {T)J 


T A 0 


: D + 1 


S (T—s)A 0 f&odW°\ 

^ odW s J ’ 


An — 


—9o 9o 
9 -9 
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Since (12) is linear, (xq(T), xn{T)) is jointly Gaussian and can be completely char¬ 
acterized by its mean and covariance matrix. We note that (—1, — 1) T is in the null 
space of Aq and thus 


E fem) = ' TA " 


-1 

-1 


-1 

-1 


In addition, Aq has the following eigen-decomposition: Aq = QqAqQq , where 


An — 


o o \ _ e (i - e f 

o -(o Q + e))' 140 00 + n 1 


i 


i Oo 

Qo 1 = Ui ! 


Then the covariance matrix is 

(52) 

1 


( Var.T 0 (T) Cav(x 0 (T),x(T)) 

\Cov(x 0 (T),x(T)) Var x{T) 

T 

o / e 


= ^Qo ^ e (T_s)A °Qg 1 ° 2 ) (QcYe {T - s)Ao dsQ 


1 


— ^Qo^Qq ) 


with 

S = 


g^gi-^ + Ooc r 2 /6»)[l-e T(e ° +e) ]\ 

_T_(_ CTo 2 + 0oa 2 m e -neo +S)] __X__ ( (jp + a 2 )[l - g-2T(0 o +9)] I . 


nai+ow/o 2 ) 


When the terminal time T is large, we can separate the middle matrix in 
the principle term and the correction term: 


into 


S = 


T(a 0 2 + 9Y/° 2 ) O' 


0 


0 


+ 


0 


^(-a 2 +0 o a 2 /0)[l -e-TVo+o)^ 

^(-a> + e 0 aye)[ 1-e-W)] + a 2 ){l ^ e -2 t (S o+6)] j ■ 


Then we have the approximation of the covariance matrix: 


(53) 

( Varxo(T) Cov(r 0 (T),x(T)) 
\Cov(x 0 {T),x(T)) Var x{T) 


1 o (nv 2 o+o 2 0 v 2 /o 2 ) 0\ T 
7V y ° l 0 n) g o 


0 


N (0 O + 0) 2 


1 1 
1 1 


From (53) we conclude that Xq(T) and x(T) are approximately equal as T becomes 
large and the probability in ( |35| is approximately P(a:o(T) £ (1,1 + da;)), which 
gives the desired rate of decay by using the fact that xq(T) is Gaussian with mean 
— 1 and approximate variance Vara^T) in (53) for large T. 


Appendix B. Proof of Proposition M 

We prove it in three steps. The first step is to show that there exists a uniform 
lower bound for J over all feasible (j). 
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Lemma 8. If h = 0, then for all <p(t,dx) such that (cf>(t,dx),x) = x{t), 

1 ( T 

j{fx 0 {t),(j){t,dx)) te[0}T ]) > —^ ( x o + h o Vo(x o )+0 o (x o -x)) 2 dt 

Za o Jo 

1 r T 

+ 2a*J ( x + 0 ( x ~ x o)) 2dt , 

for cto > 0 and for gq = 0 , 

1 f T 

j{{x 0 (t),(j)(t,dx)) te[0tT] ) > — J (x + 9(x — xo)) 2 dt, 

if x 0 + /loV'o(xo) + 9 0 (xo — x) = 0 or j((x 0 (t), </>(£, cfaO)t e [o,T]) = oo otherwise. 
Proof. By taking /( x) = x , we have 

rT {<ft - - 9-iclix - xo(t))fi\,f(x)) 2 


SU.p ^ 

Jo f(x):(<l>,(f'(x)) 2 )^0 (0i (/ (**')) ) 

f(x)=x rT (fa _ I a 2 (j) xx - 0-!^[{x - X 0 (t))(j)\,x} 2 

Jo (0i 1) 

Then we have the desired results. 


-dt 


dt = (x + 9(x — Xq )) 2 dt. 


□ 


We then prove that ff ((xo(t),p(t, da;))te[o,T]) = I\( x o(t), a:(t))g[o,T]) and conse¬ 
quently X(a;o, x) = I(xq,x). 


Lemma 9. Let p defined in |37p and h = 0. Then J({xp{t),p(t, dx))t e \p^T]) = 
I(x 0 ,x) in (23) for g 0 = 0 and j((x 0 (t),p(t,dx)) te [ 0t T]) = I{x 0 ,x) in (33) for 
g o > 0. Consequently, p(t,dx) is a minimizer and I(xq,x) = I(x q,x) for either 
gq = 0 or gq > 0 . 


Proof. By using the same argument in ma Proposition 5.3], if <j)(t,dx) is absolutely 
continuous with respect to the Lebesgue measure with the smooth density function 
x), then 


sup 

fC):(<l>,(f'(x)) 2 )^ o 


{(ft - 2 G 2 (f xx - 0fh[{x - x 0 (t))<j>\,f{x)y 




dt = 


(</>, (g(t,x)) 2 )dt, 


where g(t, x) satishes 

\ d d 

4>t(t,x) - -G 2 (f xx (t,x) - 9— [(a; - x 0 (t))<f>(t,x)\ = — {<f(t,x)g{t,x)). 

li<p{t,x) = p(t, x), then by using the fact that pt = —x(t)p x and \G 2 p xx +9-j((^[{x — 
x o(t))p ] = @[(x(t) — x o{t))p x ], the corresponding g(t,x) satisfies 

, _ _ _ d 

x{t)p x - 9[(x(t) - Xq( t))p x \ = —(g(t,x)p). 

Then g(t, x) = —x(t) — 9(x(t) — x 0 (t)) and / O T (0, (g(t, x)) 2 )dt = f X (x(t) + 9(x(t) — 
Xq( t))) 2 dt. We therefore obtain the desired results. □ 


Finally we show that the minimizer (p(f, dx)) tg [o,T] is unique. 

Lemma 10. The minimizer (p{t, dx)) te [ 0tT ] of inf^p^) j({x 0 {t), dx)) te [ 0>T ]) 

is unique for all {<j>{t, da;))te[o,T] such that {</>(t,dx),x) = x{t) for all t £ [ 0 ,T] and 
0 ( 0 , dx) = p( 0 , dx). 
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Proof. From the previous lemmas we conclude that if {<f(t, daO)te[o,T] is a minimizer, 
then 

{(ft - \° 2 (fxx - x 0 (t))cj)}J{x )) 2 

x= argsup --- , , , -• 

Therefore for any perturbation f(x), 

{(ft - \<J 2 <t>xx - S-^[(x - x 0 (t))(/)\,x + ef(x)} 2 _ o 

(<f>,(l + ef'(x)y) 

which leads to 

\ () \ Q 

(ft - 2 <j2( f >xx ~ d ~gf^ x ~ x ° W)^ 1 ] = (*At ~ 2 cr2( t >xx ~ ^dx^ X ~ x ° (*))$]> 

= [x(t) + 9(x(t) - x 0 (t))]<f. 

In other words, a minimizer {(f{t,dx)) te i 0 ,t] must satisfy the above linear parabolic 
PDE that has a unique solution with the given initial condition </>(0, dx) = p{ 0, dx). 

□ 


d_ 

de 


Appendix C. Proof of Proposition [7] 


We can rewrite the problem in the matrix form: 


min -E 

( a (*))te[o,T] 2 


f a(t) T Ra(t) + X(t) T QX(t)dt 
Jo 


dX = 'ZdW + AX + Badt, 


where 


' _£o_ n,.T 

S = ( VN Um 
0m crl 


Q = 9 C 


N -« T 
—u I 


A = 

, R = 


-0 o -H o 

Ou -01 

1 / 1 Oit T 


B = 


0 0m t 
Ou I 


u= (1,...,1) T . 


9 C \0u I 

We apply the standard theory [23] Theorem 6.1] and we find that the optimal 
control is 

aft) = -Rr^S^X^) 

where S ft) is solution of the matrix Riccati equation 


- — S = A t S + SA - S t BR~ 1 B t S + Q, 

at 

with the terminal condition S(T) = 0. We find that 

cu\-f Na ^ b ^ uT \ 

\ bff)u d(f)I+^J/’ 


where J is the N x N matrix full of ones and (a(f), 6(t), d(t), e(f)) tg [ 0> T] is the 
solution of 


aft) = 2(0 O + Ho)a(t) — 2 9b(t) + 9 c b 2 ft) — 9 C , 

bft) = (0o + H 0 + 9)b(t) — 9d(t) — 0 o a(f) + 9 c bft)dft) + 0 C — 9 eft) + 9 c bft)e(t), 

dft) = 29d(t) + O c d 2 (t) — 0 C , 

eft) = —2 9 0 b{t) + 2 9e(t) + 9 c (2dft)e(t) + e 2 (t)), 
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with (a(T), b(T),d(T), e(T)) = (0, 0,0, 0). Therefore the optimal control is 
= -0 c (b(t)X Q (t) +d(t)Xj(t) + e(t)X N (t)), j = 
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